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Abstract 

Modelling stochastic systems has many important applications. 
Normal form coordinate transforms are a powerful way to untangle 
interesting long term macroscale dynamics from detailed microscale 
dynamics. We explore such coordinate transforms of stochastic differ- 
ential systems when the dynamics has both slow modes and quickly 
decaying modes. The thrust is to derive normal forms useful for macro- 
scopic modelling of complex stochastic microscopic systems. Thus we 
not only must reduce the dimensionality of the dynamics, but also 
endeavour to separate all slow processes from all fast time processes, 
both deterministic and stochastic. Quadratic stochastic effects in the 
fast modes contribute to the drift of the important slow modes. The 
results will help us accurately model, interpret and simulate multiscale 
stochastic systems. 
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1 Introduction 



Normal form coordinate transform ations provide a sound basis for simplify- 
ing multiscale nonlinear dynamics (jElphick et al.lll987l . ICox fc Robertslll995 . 
e.g.). In systems with fast and slow dynamics, a coordinate transform is 
sought that decouples the slow from the fast. The decoupled slow modes 
then provide accurate predictions for the long term dynamics. Arguably, 
such normal form coordinate transformations that decouple slow and fast 
modes provide a much more insightful view of simplifying dynamics than 
other, more popular, techniques. A veraging is per haps the most popular 



technique for simplifying dynamics (jVerhulstl 12005k Ch apters 11-13, e.g.), 



espec i ally for stochastic dynamics that we explore here (IPavliotis fc Stuart 



20061 . iGivon et al.l 120061 . e.g.). But averaging fails in many cases. For ex- 



ample, consider the simple, linear, slow-fast, system of stochastic differential 
equations (sdes) 



dx = ey dt and dy = — y dt + dW , 



(1) 



where for small parameter e the variable x(t) evolves slowly compared to 
the fast variable y[t). Let us compare the predictions of averaging and a 
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'normal form' coordinate transform. First consider averaging: the fast vari- 
able y, being an Ornstein-Uhlenbeck process, rapidly approaches its limiting 
PDF that is symmetric in y. Then averaging the x equation leads to the 
prediction dx = dt ; that is, averaging predicts nothing happens. Yet the 
slow x variable must fluctuate through its forcing by the fast y. Second, and 
similar to illuminating coordinate transforms used in Sections [2HH1 introduce 
new coordinates X and Y to replace x and y where 



X- eY 



e~ x dW T and y = Y + 



_t dW T 



In the X and Y coordinates the SDE system ([T]) decouples to simply 



dX = edW and dY 



-Y dt 



(2) 



(3) 



In these new coordinates Y — > exponentially fast. Thus in the long term 
the only significant dynamics occurs in the new slow variable X which sys- 
tem shows undergoes a random walk. The method of averaging com- 
pletely misses this random walk: true, the mean x remains at zero; but the 
growing spread about the mean is missed by averaging. Stochastic coordi- 
nate transforms such as (j2J) decouple fast and slow variables to empower us to 
extract accurate models for the new slow variable X. They are called 'normal 
form' transformations because this decoupling of stochastic dynami cs is anal- 



ogou s to correspon ding simplifications in deterministic systems ( iMurdock 
20031 . lArnoldl 120031 . e.g.). This article establishes useful properties for such 



stochastic normal form coordinate transformations in modelling multiscale 
nonlinear stochastic dynamical systems. 

Stochastic odes and pdes have many important applications. Here we 
restrict attention to nonlinear SDEs when the dyna mics of the SDE has both 
long lasting slow modes and decaying fast modes (lArnold fc Imkellerlll998 . 
e.g.). The aim underlying all the exploration in this article is to derive normal 
forms useful for macroscopic modelling of stochastic systems when the sys- 
tems are specified at a detailed microscopic level. Th us we endeavour to sep- 
arate all fast t ime processes from the slow processes ( IChao &: Robertslll996 . 
Robert si 12006 a e.g.). Such separation is especially interesting in stochastic 
systems as white noise has fluctuations on all time scales. In contrast, al- 
most all previous approaches have been content to derive normal forms that 
support reducing the dimensionality of the dynamics. Here we go further 
than other researchers and additionally separate fast time processes from the 
slow modes. 



Arnold fc Imkellerl (119981 ) developed rigorous theory for stochastic coor- 
dinate transforms to a normal form in the "non-resonant" cases. Then a 
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stochastic system is rigorously equivalent to a decoupled linearised system 
no matter how disparate the time scales. However, sto chastic fast-slow sys- 
tems are essentially resonant. lArnold fc Imkellerl (119981 ) report, §4, that their 
results do app ly, but that the resulting normal form is generically nonlinear 
(lArnoldl 120031 . §8.4 also). They comment that the normal form transfor- 
mation involves anticipating the noise processes, that is, involving integrals 
of the noise over a fast time scale of the future. However, in contrast to 



the examples of I Arnold fc Imkellerl (119981 ) (jArnoldll2003l . corrected), Sections 
[2] and [3] argue that such anticipation can always be removed from the slow 
modes with the result that no anticipation is required after the fast transients 
decay. Furthermore, Sections [2] and [3] argue that on the stochastic slow man- 
ifold all noise integrals can be removed from terms linear in the noise to leave 
a slow mode system, such as the simple dX = e dW of the normal form Q, 
in which there are no fast time integrals at all. The arguments demonstrate 
that, except for some effects nonlinear in the noise, all fast time processes 
can be removed f rom the slow modes of a n ormal form of stochastic systems. 

The theory of lArnold fc Imkellerl (11998 ) appli es only to finite dimensional 
stochastic systems. Similarly, |Pu fc Duanl (|2006|)'s theory of invariant man- 
ifold reduction for stochastic dynamical systems also only applies in finite 
dimensions. But many applications are infinite dimensional; for example, the 
discretisation of stochastic pd es approximates an inertial manifold of stochas- 
tic dynamics (lRobertsll2006d). F ollowi ng the wide recognition of the u tility of 



inertial manifolds (]Temamlll990l . e.g.). iBensoussan fc Flandolil (119951 ) proved 



the existence of attractive stochastic inertial manifolds in Hilbert spaces. The 
stochastic slow manifolds obtained in Sections [2] and [3] via stochastic normal 
forms are examples of such stochastic inertial manifolds, albeit still in finite 
dimensions. 

To derive a normal form we have to implement a coordinate transfor- 
mation that simplifies a stochastic system. But the term 'simplify' means 
different things to different people depending upon how they wish to use 
the 'simplified' stochastic system. Our aim throughout this article is to cre- 
ate stochastic models that may efficiently simulate the long term dynamics 
of multiscale stochastic systems. This aim is a little different to that of 
previo us researchers and so the resu lts herein are a litt l e diffe rent. For ex- 
ample, ICoullet &: Spiegell (119831 ) and lArnold &; Imkellerl (119981 ) do not avoid 
fast time integrals because their aim is different. Principles that we require 
are the following: 

1. Avoid unbounded (secular) terms in the transformation and the evolu- 
tion (ensures uniform asymptotic approximations); 
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2. 



Decouple all the slow processes from the fast processes (ensures a valid 
long term model); 



3. Insist that the stochastic slow manifold is precisely the transformed 
fast modes being zero; 

4. Ruthlessly eliminate as many as possible of the terms in the evolution 
(to simplify at least the algebraic form of the SDEs); 

5. Avoid as far as possible fast time memory integrals in the evolution (to 
endeavour to remove all fast time processes from the slow modes). 

In general we can meet all these principles, although the last two are easy 
as they are only phrased as 'as far/many as possible': Section [2] explores the 
issues in a particular example stochastic system; whereas Section [3] presents 
general arguments for finite dimensional, nonlinear, stochastic differential 
systems. Sections 12 . 31 and [2.41 additional show that two other alternative, and 
superficially attractive, principles are not useful in the context of macroscale 

mo delling of stochastic systems . 

Sri Namachchivaya fc Leng ( 199dl ) and Sri Namachchivaya fc Lin ( 1991 ) 



emphasise the importance of effects quadratic in the stochastic noise "in 
order to capture the stochastic contributions of the stable modes to the drift 
terms of the critical modes." Sections [2] and [3] also address such important 
quadratic effects. The generic result of this normal form approach is that not 
all the memory integrals can be removed from the evolution of the stochastic 
slow variables: some terms quadratic in the noise retain memory integrals. 

Section H] explores the implications of these results for macroscale sim- 
ulation of stochastic systems. The normal form approach empowers us to 
address the effect of anticipatory integrals, the influence of the noise on av- 
erages, especially noise induced drift, and the failure of averaging to provide 
a systematic basis for macroscale simulation. 

Lastly, Section [3] discusses in detail a normal form of a stochastically 
forced Hopf bifurcation. Not because it is a Hopf bifurcation, but instead 
because it is a generic example of stochastic effects in oscillatory dynamics. 
The primary issue is how to 'average' over both the nonlinear oscillation and 
the noise effects to generate a prescriptive model of the dynamics over much 
longer time scales. A complex valued, time dependent, coordinate transform 
can, with considerable care, derives a model SDE that is valid for simulat- 
ing the long term evolution of the stochastic oscillating dynamics. A future 
application coul d be to the m odelling of atmospheric white noise forcing of 



oceanic modes: iPiercd ( 120011 ) discusses this situation from an oceanogra- 
pher's perspective. 
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Figure 1: sample trajectories of the example stochastic system (jSJ) from dif- 
ferent initial conditions for different realisations of the noise, a = 0.05 . 



2 Explore in detail a simple nonlinear stochas- 
tic system 

This section considers the dynamics of one of the most elementary, nonlinear, 
multiscale stochastic systems: 



dx = — xy dt and dy 



-y + x 2 -2y 2 )dt + adW. 



(4) 



The issues raised, and their resolution, in this relatively simple stochastic 
system are generic as seen in Section [31 

Thr oughout this article I ado pt the Stratonovich interpretation of SDEs, 
as does Arnold fc Imkellerl <jl998h . so that the usual rules of calculus apply. 
To ease asymptotic analysis I also adopt hereafter the notation of applied 
physicists and engineers. Thus I formally explore the SDE system (j4|) in the 
equivalent form of 



x 



-xy and y = — y + x 2 — 2y 2 + a4>(t) , 



(5) 



where overdots denote formal time derivatives and the 'white noise' cb(t) is 
the formal time derivative of the Wiener process W(t) . Both the Stratonovich 
interpretation and the adoption of this formal notation empowers the use of 
computer algebra to handle the multitude of details. 
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Figure [T] plots some typical trajectories of the SDE system (J3J). In this 
domain near the origin the y variable decays exponentially quickly to y m 
x z ; whereas the x variable evolves relatively slowly over long times. Thus 
the y variabl e represents fast, micro scopic, uninteresting modes, they are 



"slaved" to x (jSchoner &: Hakenlll986l ). whereas the x variable represents the 



long lasting, macroscopic modes of interest to the long term dynamics. The 
white noise crc|>(t) added to the y dynamics induces fluctuations on all time 
scales. We explore fundamental issues in the modelling of the multiscale 
dynamics of the SDE (j^J). 

The challenge is to adapt the deterministic normal form transformation, 
Section 12.11 to the stochastic system (jSJ) in order to not only decouple the 
interesting slow modes, but to simplify them as far as possible, Section 12.21 
Two other alternative normal forms are explored in Sections 12.31 and 12.41 
but I argue that they are not so useful. The analysis and argument is very 
detailed in order to demonstrate in a simple setting how Principles [TH5] are 
realised at the expense of having to anticipate future noise. If you are familiar 
with the concept of stochastic normal forms, you could skip to Section [3] for 
generic arguments of the new results. 

2.1 Decouple the deterministic dynamics 

Initially consider the example toy system ([5|) when there is no noise, cr = . A 
deterministic normal form coordinate transform decouples the deterministic 
slow and fast dynamics: 

x = X + XY + |XY 2 - 2X 3 Y + f XY 3 + • • • , (6) 
y = y + X 2 + 2Y 2 + 4Y 3 - 4X 2 Y 2 + 8Y 4 + ■ • • . (7) 

Figure [2] shows the coordinate curves of this (X, Y) coordinate system. The 
coordinate transform is a near identity because near the origin x « X and y « 
Y. In the new (X, Y) coordinate system, the evolution of the toy system (jSJ) 
becomes 

X = -X 3 and Y = -(1 + 2X 2 + 4X 4 )Y + • • • . (8) 

Observe the Y-dynamics are that of exponentially quick decay to the slow 

manifold Y = at a rate 1 + 2X 2 + 4X 4 H . In the original variables, 

from th e transform (El) and ((Jj), this slow manifold is the curve x = X and 



y = X 2 (lElphick et al.lll987l ). The dynamics on this slow manifold, X = — X 3 
from (JBD, form the accurate, macroscopic, long term model. 

The slow X dynamics are also independent of the Y variable and thus 
the initial value Y(0) and subsequent Y(t) are immaterial to the long term 
evolution. Thus to make accurate forecasts, project onto the slow manifold 
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Figure 2: coordinate curves in the xy-plane of the (X,Y) coordinate sys- 
tem that simplifies to (jSJ) the algebraic description of the dynamics of the 
deterministic (cr = 0) system (j^J). 



S 



Y = along the coordinate curves of constant X flCox fc Robertslll995l ) seen 
in Figure [21 Equivalently, because the slow X dynamics are independent of 
the Y variable, the dynamics of the system (jSJ) map the curves of constant X 
in Figure [2] into other curves of constant X. Thus initial conditions on any 
one curve of constant X all evolve towards the same trajectory on the slow 
manifold. 

But these comments are all for deterministic dynamics, ct = . The 
next subsection answers the question: how can we adapt this beautifully 
simplifying coordinate transform to cater for stochastic dynamics 



? 



2.2 Simplify stochastic evolution as far as possible 

Now explore the construction of a coordinate transform that decouples the 
fast and slow dynamics of the toy SDE (JHJ) in the presence of its stochastic 
forcing. In order to cater for the stochastic fluctuations, the coordinate 
transform must be time dependent through dependence upon the realisation 
of the noise, as shown schematically in Figure [3J This subsection is very 
detailed in order to argue that no alternatives go unrecognised. The method 
is to iteratively refine the stochastic coordinate transform based upon the 
residuals of the governing toy sde flSJ). 

Although our focus is on the case when cb(t) is a white noise, because we 
use the usual calculus of the Stratonovich interpretation, the algebraic results 
also apply to smoother processes cb(t). For two exa mples, the forcing 4>(t) 



could be the output of a deterministic chaotic system (lJust et al.ll2001l . e.g.), 
or the forcing cj)(t) could be even as regular as a periodic oscillator. Thus the 
algebraic expressions derived herein apply much more generally than to just 
white noise cj). However, the justification for the particular coordinate trans- 
form often depends upon the peculiar characteristics of white noise. For 
forcing cj) which is smoother than white noise, although our results herein 
apply, other particular coordinate transforms may be preferable in order to 
achieve other desirable outcomes in the transformation (outcomes not attain- 
able when cj) is white noise). These possibilities are not explored. Instead, 
almost everywhere throughout this article, the forcing cb(t) denotes a white 
noise process in a Stratonovich interpretation of SDEs. 

Let us proceed to iteratively deve lop a stochasti c coordinate transform of 



the SDE (jHJ) via stepwise refinement ([Roberts! 119971 ). 



First, consider the fast dynamics With x rs X and X « 0, seek a 
change to the y coordinate of the form 

y =Y + ri / (t,X,Y) + --- and Y = — Y + G'(t, X, Y) + ■ • • , (9) 



9 



i — i — i — i — i — i — i — i — i — i — i — i — i — i — i — i — r 



-0.5 -0.4-0.3-0.2 -0.1 0.0 0.1 0.2 0.3 0.4 0.5 

X 



Figure 3: four different (coloured) meshes represent either four realisations 
sampled at one instant or one realisation sampled at four instants of the 
coordinate curves in the xxj-plane of the stochastic (X, Y) coordinate system 
that simplifies to 02 6 j) and (J29l) the algebraic description of the dynamics of 
the stochastic (a = 0.2) system ([5]). 
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where r\' and G' are small, 0(e 2 ), corrections to the transform and the 
corresponding evolution. I introduce the parameter e to provide a convenient 
ordering of the terms that arise in the algebra: formally set e = |(X,Y, <r)| 
with the effect that e counts the number of X, Y and a factors in any one 
term. Substitute (jHD into the y SDE (JHJ) and drop products of small corrections 
to recognise we need to solve 

G / + ^+t 1 '-Y^ = ac( ) (t)+X 2 -2Y 2 +.--; (10) 

partial derivatives are here done keeping constant the other two variables of 
X, Y and t. 

Consider first the deterministic terms. To solve G' + r\' t +V — Yt)y = 
X 2 — 2Y 2 , keep the evolution as simple as possible (Principle H|) by seeking 
corrections G' = and r\' = aX 2 + bY 2 . Substitute into the equation to see 
a = 1 and b = 2 . Consequently include r\' — X 2 + 2Y 2 into the coordinate 
transform. 

Second, consider the stochastic term 0"<b(t) in the right-hand side of (11111) . 
To solve G' + r\' t + r\' — Yqy = 0"4>(t) keeping the Y dynamics as simple as 
possible (Principle Hj) choose the convolution cre^^cf) , defined in ( I12j) . to be 
part of the correction r\' to the coordinate transform. Consequently the new 
approximation of the coordinate transform and the dynamics is 

y =y + X 2 + 2Y 2 + a-e~ t ★$ + ■■• and Y = — Y + ■ ■ • . (11) 

In these leading order terms in the coordinate transform, see the stochastic 
slow manifold (ssm) Y = corresponds to the vertically fluctuating parabola 
y ~ X +ae _t ^4' as seen m the overall vertical displacements of the coordinate 
mesh in Figure [3j 



The convolution For any non-zero parameter and consistent with the 
convolution in the example transform (j2J), define the convolution 

j J^exp^t-Tfl^MdT, [i<0, 
^expMt-T)]*(T)dT, [i>0, 1 ; 

so that the convolution is always with a bounded exponential (Principle [T]). 
Five useful properties of this convolution are 

e-*1=ilj. (13) 

^-e^^cj) = — sgn(x4> + [le^-kfy , (14) 
dt 
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E[e^*c|>] = e^*E[<|>], (15) 
E[(e^*cM 2 ] = -L (16) 



l 



e Ht * e vt * = <^ L ' J ' r " ' ~ " ' (17) 

|_2Eiii[e^*-e vt *] , ^iv>0&|x^v. v ; 

Also remember that although with \x < the convolution e^ 1 * integrates over 
the past, with \x > 0, as we will soon need, the convolution e Ht * integrates 
into the future; both integrate over a time scale of order 1 



Second, consider the slow dynamics Seek a correction to the stochastic 
coordinate transform of the form 

x = X + £/(t,X,Y) + 0(e 3 ) and X = F'(t, X, Y) + £>(e 3 ) . (18) 

where £,' and F' are C*(e 2 ) corrections to the transform and the corresponding 
evolution. Substitute into the x equation of SDE (jSJ) and omit small products 
to recognise we need to solve 

F' + ^ - = -XY + o-Xe- 1 ★ 4> + G(e 3 ) . (19) 
ot or 

First, try £,' = aXY to find the deterministic term on the right-hand side is 
matched when a = 1 . Second, consider the stochastic part of the equation: 
F' + £,{ — Y£,y = crXe _t * (J) . The E,y cannot help us solve this stochastic part 
as there is no Y factor in the right-hand side term. We do not want to assign 
a fast time convolution into the slow evolution F' (Principle [5]), but cannot 
integrate the forcing <b into £,' as then terms would grow like the Wiener 
process W = J cp dt (Principle [1]). Instead, formally integrate by parts to 
split e _t * <J) = — 4> + e _t * <b and hence choose components F' = — crXcp and 
£/ = crXe _t * (b . Consequently 

x = X + XY+aXe" t *ct) + C(e 3 ) and X = — crXcp + C(e 3 ) . (20) 

Third, reconsider the fast dynamics Seek corrections, r\' and G', to the 

y transform and Y evolution driven by the updated residual of the y equation 
in SDE flSJ: 

G' + ^+ti'-Y^ = -4aYe- t ^4)-2o- 2 (e- t ^4)) 2 + C(e 3 ) . (21) 
Separately consider the two stochastic forcing terms on the right-hand side. 
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• To solve G' + r\' t + r\' — Yr\y = — 4crYe~ t * cj) we must seek G' and r\' 
proportional to Y, whence r\' — Yr\' = , and so to avoid secular terms 
in r\' (Principled]), and to avoid fast time convolution in the Y evolution 
(Principle EJ), integration by parts enables to choose G' = — 4crY(|) and 
rj' = 4o-Ye _t * 4> . 

• To solve G'+ri^+ri' — Yriy = — 2cr 2 (e~ t *cb) 2 seek G' andr)' independent 
of X and Y. Hence choose G' = (Principle H|) and then the convolution 
T| ' = — 2cr 2 e~ t * (e _t * <b) 2 corrects the coordinate transform^] 

Consequently, the fast time transform and dynamics are more accurately 

y = Y + X 2 + 2Y 2 + a[e- t *4)+4Ye- t *cb] 

-2a 2 e- t *(e- t *c[)) 2 + C(e 3 ) , (22) 
Y = -Y-4aYct) + C(e 3 ) . (23) 

Fourth, reconsider the slow dynamics Seek corrections to the trans- 
form and evolution, £/ and F, driven by the updated residual of the x equa- 
tion of SDE 

F' + ^ - Y^y = -X 3 - 3XY 2 + aXY(54) - Ge^ * 4>) 

+ a 2 X [4>e _ * * 4> — (e^* cp) 2 + 2e _t *(e _t *cb) 2 ] +C(e 4 ). (24) 

Consider the right-hand side term by term: 

• To account for the deterministic forcing, we must choose F' = —X 3 and 
£/ = |XY 2 in the traditional manner. 

• To match the term linear in noise, F' + L,' t — Y£,y = crXY(5(|) — 6e~ t *cb) , 
we must seek F' and £,' proportional to XY, whence L,[ — YE,' y M L,[ — 
£/ . Consequently foreknowledge, anticipation, of the noise appears. 
Consider the two cases: 

— allow ing anticipation (implem enting Principle Hj) and in accord 



with lArnold Sz Imkellerl (119 98). we assign all of this term to the 
coordinate transformation with F' = and £,' = — o"XYe +t * (54> — 
6e _t -k 4>) ; 



lr The right-hand side of this correction r\' has non-zero mean. We could assign the 
mean, — cr 2 , into the Y evolution as a mean (downwards) forcing, but then this destroys 
Y = as being the slow manifold, contradicting Principle [3] 



13 



— disallowing anticipation, we must assign all of this term into the 
X evolution by assigning F = crXY(54> — 6e~ t * cb) and £/ = 
— the difficulty here being that the evolution to the SSM then de- 
pends undesirably upon Y, contradicting Principle El Section 12.41 
explores this case. 

• For the quadratic noise term in (I24p . seek contributions to the solution 
which are proportional to X; consequently, on the left-hand side — Y£,y = 
. Then we avoid secularity (Principle [1]) by extracting the mean of the 
right-hand side term and assign the mean into the evolution F; but at 
least part of the fluctuations cannot be assigned into the transform £,' as 
the integral of noise is a Wiener process which almost surely is secular. 

Now, as in the earlier integration by parts, separate these quadratic 
noise terms into 

(e-^ct.) 2 = cbe-Hcb-^tle-^cb) 2 ], 
e^*(e^*<|>) 2 = (e-^cb) 2 -^-**^ - ***) 2 ] 

= cbe^ * cb - f^le- 1 * cb) 2 + e"* * (e"* * cb) 2 ] , 

and so these contribute corrections F = o^Xcbe -1 *^ and £,' = — o" 2 X(^+ 

The upshot is that the x transformation and X evolution is more accurately 

x = x + XY+fXY 2 + cr[Xe- t *ct. + XYe +t *(-54) + 6e- t *4))] 

-a 2 X(I+2e-V)(e- t *c))) 2 + C(e 4 ) , (25) 
X = -X 3 -ffXcb + 2cr 2 Xcbe- t *cb + e>(e 4 ) . (26) 

Lastly, reconsider the fast dynamics Seek corrections r\' and G' to the 

y transform and Y evolution driven by the updated residual of the y equation 
in SDE (ED: 

G' + ^+V-Y^ = -2X 2 Y-8Y 3 
ot or 

+ cr[2X 2 (cb - e _t * cb) + 8Y 2 (2cb - 3e _t * cb)] 

+ S^Ypcbe - * * cb - l[z~ x * cb) 2 + e _t * (e _t * cb) 2 ] 

+ 8ff 3 (e~ t *cb)e- t *(e~ t *cb) 2 . (27) 

Separately consider the forcing terms on the right-hand side. 

• The deterministic terms generate contributions G' = — 2X 2 Y and r\' = 
4Y 3 in the usual manner. 
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• To solve G' + T|( + X[' — Yt|y = +2aX 2 (4) — e _t * <b) we must seek 
G' and r\' proportional to X 2 whence — Yr|' = 0. Consequently, assign 
everything to the transform employing the convolution e _t *: G' = 
and T] ' = 2crX 2 (e~ t * (J) — e _t * e _t * cb) (Principle [3] and Principle HJ) . 

• The Y 2 term raises the issue of anticipation again. To solve G' +T]{ + 
r\' — Yt|y = 8crY 2 (2cb — 3e _t * cb) we must seek solutions proportional 
to Y 2 and hence the left-hand side G' + r|{ + r\' — Yn Y = G' +r\' t — r\' . 
By Principle HI assign the entire forcing into the transform by setting 
G' = andr]' = — 8<rY 2 e +t *(2cb — 3e _t *c()) . This requires anticipation 
of the forcing through e +t * (b and e +t * e _t * cb . 

• The remaining two terms do not generate new issues so I do not describe 
the details. 

The normal form coordinate transform for the fast dynamics is thus 

y = Y + X 2 + 2Y 2 + 4Y 3 + cr[e- t *ct. + 4Ye- t *ctj 

+ 2X 2 e~ t * (cb - e _t * cb) - 8Y 2 e +t * (2(f) - 3e _t * cb)] 
+ 4cr 2 Y(1 -2e^*)(e^*<|>) 2 

+ 8o-V t * [(e^cbje-Me^cb) 2 ] +C(e 4 ) , (28) 
Y = -Y-2X 2 Y-4o-Yc|) + 8o- 2 Yct)e" t ^c|) + C(e 4 ) . (29) 

Higher order model Further algebra leads to the construction of a stochas- 
tic coordinate transform from the original [x,y] variables to the new (X, Y) vari- 
ables so that the dynamics of the example SDE (JHJ) is 

X = -X 3 - crXcp + 2a 2 X4)e- t * 4> 

- 4o- 2 X 3 4)e^ * e _t * cb + O (e 6 , a 3 ) , (30) 

Y = -(1 +2X 2 + 4X 4 )Y-4a(1 + X 2 )Ycb + 8cr 2 Ycpe * * cb 
+ 4cr 2 X 2 Ycb [3e _t * cb - e +t * 4) - 2e _t * e _t * cb] 
+ 0(e 6 ,cr 3 ). (31) 

By employing a coordinate transform that depends upon the realisation of 
the noise, we maintain Y = as the exponentially attractive SSM, see ()3ip . 
independent of realisation!! However, we can only make the X evolution 
independent of Y, as in flHUj) . by anticipating noise, albeit anticipating only 
on a fast time scale into the future. To rationally project onto the SSM we 
must accept some fast time scale anticipation. 

2 Possibly, large enough noise cr could make Y = repulsive, just as large enough X may 

do. 
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Irreducible fast time convolutions generate drift Also note that the 
X and Y evolution equations, (13"U1) and f[3"T]) . contain algebraically irreducible 
nonlinear noise such as <pe^ 1 * <j), in defiance of Principle [5j Over long 
times such irreducible noise could be replaced by \ + -^<J) for some effec- 
tively new noise <J)(t) f Chao fc Roberts! Il996l ). Such replacement was also 
j ustifi ed by Khasminskii (1996) as described by ISri Namachchivaya fc Leng 
(jl990l ). Importantly, such quadratic noise, in effect, generates a mean deter- 
ministic drift term in the slow dynamics ( Sri Namachchivaya fc Leng 199dl . 
Sri Namachchivaya fc Linlll99ll . e.g.). In applications such drifts can be vital. 



The average SSM is not the deterministic slow manifold For the 

toy sde Section |2~T1 shows the deterministic slow manifold is y = x 2 . In 
general the SSM fluctuates about a mean location which is different to this 
deterministic slow manifold. From ( 1251) and (}28l) with fast variable Y = , 
the SSM is 

x = X + o-Xe- t *4)-o- 2 X(^ + 2e-V)(e- t ^ct)) 2 + C(e 4 ) , (32) 
y = X 2 + o-[e~ t *c|) + 2X 2 e- t *(*-e~ t *cb)] 

+ 8o 3 e- t * [(e-^c^e-Me^ct)) 2 ] +C(e 3 ) , (33) 

Take expectations, and using (1151) and (1161) . 

E[x] = (1 - fcr 2 )X + C(e 4 ) and E[y] = X 2 + O(o- 3 , e 4 ) . (34) 

Observe E[u] ^ E[x] 2 , instead Ehj] « (1 + |cr 2 )E[x] 2 so that the average 
SSM is a steeper parabola shape than the deterministic slow manifold. It is 
quadratic noise processes that deform the average SSM from the deterministic. 

Forecast from initial conditions Suppose at time t = we observe the 
state [xo,yo), what forecast can we make with the SSM SDE ( 130]) ? Revert the 
asymptotic expansion of the stochastic coordinate transform ( 1251) and ( 1281) 
to deduce 

X = x + x 3 — xy + |xy 2 + 2oxue +t * <j) — 2cr 2 x(e +t * <|))(e _t * 4>) 

+ 0{e 4 ), (35) 
Y = n-x 2 -2u 2 -ae- t *4) + (T 2 (l +e- t ^)(e- t *4)) 2 + C(e 3 ) .(36) 

Then the correct initial condition for the long term dynamics on the SSM, 
governed by the SDE ( JHUj) . is the X component of this reversion, 035p . evalu- 
ated at the observed state, namely 

X(0) = x +Xo-xoUo+|xono+2CTx yoe + Vc()-2o- 2 xo(e + V4))(e _ Vct))+C'(e 4 ) . 

(37) 
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This is a projection of the observed initial state onto the SSM to provide an 
initial condition X(0) for the slow mode. However, this projection involves 
both memory and anticipatory convolutions of the noise. There are at least 
three interesting issues with computing this initial X(0). First, at the initial 
instant we do not know either the future nor the past, so the terms involv- 
ing the noise cb are unknown. Using the expectations (1151) and (FIBj) . the 
projection X(0) has known mean 

E[X(0)] = x + x^ - x oVo + \xovl + ^(e 4 ) , 

with known variance 

V[X(0)]^2a 2 x 2 y 2 + a 4 x 2 . 

That is, a given observed state (xo,yo) corresponds to a stochastic state 
for the evolution of the slow mode model on the SSM. Second, but if this 
state X(0) for the slow mode is to be used in a simulation to make forecasts 
of the future, then we know the future of the noise cb. The future values of 
noise cb are just those we use in integrating the slow mode SDE ( l30i) . Thus 
for simulation, we do eventually know the anticipatory convolutions e +t * cb 
in (|37|) . but not the memory convolution e _t * cb. In this case the mean of 
the projection 

E[X(0)] =x + x^-x yo + fxoyo + 2o-x yoe +t *ct) + 0(e 4 ) , 

with variance V[X(0)] ~ 2cr 4 Xo(e +t * cb) 2 . Lastly, if we made additional ob- 
servations for times t < , then the additional information could partially 
determine the past history of the noise cb and hence help us estimate the 
memory convolution e _t * cb. These three cases emphasise that the initial 
state X(0) of the slow variable depends upon more than just the initial ob- 
served state (xo,Do)- 

2.3 Try retaining some noise forcing of the decaying 
modes 

What if, contradicting both Principle H] and Principle El we allow some forcing 
noise to remain in the Y evolution? Is there any useful freedom? Here I argue 
there is not. 

First, consider the y dynamics Assume that we solve equation (TT0|) by 
choosing r\' = and G' = crcb . This keeps the noise imposed upon the rapid 
decay of the fast variable Y. 
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Second, consider the x dynamics The residual of the x equation in 
SDE (jSj) then leads to solving 

F + £,( - = -crXct, + • • ■ , (38) 

instead of (fT9"l) . The solution must be proportional to X, so the Y£,y term on 
the left-hand side is no use. But secularity, Principle HJ implies we cannot 
assign any of the right-hand side component into £, and so we are forced to 
have £,' = and F' = — crXc() . So far the coordinate transform itself remains 
identical to the deterministic 

Third, reconsider the y dynamics Using the above coordinate trans- 
form, determine more corrections from the new residual by solving 

G'+ri( + n , -Yn; = -4o-Yct> + 2o-(X 2 -6Y 2 )c|) + -- - . (39) 

The terms on the right-hand side could all be assigned to the Y evolution as 
G' = cr(— 4Y + 2X 2 + 4Y 2 )cJ) and r\' = . The term linear in Y certainly has 
to go into the Y evolution, but the others do not. Instead, try an arbitrary 
convex combination of the above and G' = — 4crYcJ) + 8ff 2 Ycbe +t * cj) (an 
anticipatory convolution appears) and £,' = cr(2X 2 e _t * cf) — 4Y 2 e +t * cb) . 

Last, reconsider the x dynamics Seek corrections £/ and F' to the 

x transform and the X evolution driven by the new residual: 

F' + - Y£,^ = 2oXYct, + ■■■ . (40) 

This residual is independent of the convex combination of the terms just de- 
termined above by the y dynamics. The intractable difficulty here is that the 
term 2crXY4> on the right-hand side has to be assigned into F', contradicting 
Principle El Further, we also get no hint of the quadratic stochastic mean 
drift effect. Since I parametrised the freedom in the previous step, and the 
second step was forced, then to obtain the quadratic mean drift term we must 
change the first step. That is, we can only satisfy Principle [2] and also can 
only extract the quadratic mean drift term by abandoning the assumption 
that the first step is useful. Thus, allowing nonhomogeneous forcing of the 
Y evolution is not useful for the SDE (JSJ). 

2.4 Avoiding anticipation is less useful 

Alternatively, suppose we disallow anticipatory convolutions. At least one 
of the Principles [EH5] then has to be abandoned. Principles H] and [5] are "as 
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possible" principles so we meet them as best we can. Well ordered asymptotic 
expansions are essential, so avoid secularity, Principle [TJ Section 12.31 shows 
abandoning Principle [3] is ineffective. Consequently, in this subsection we 
explore abandoning Principle [2, the requirement to completely decouple the 
slow modes from the fast modes. But abandoning this principle means we 
are no longer able to use the slow model to make high accuracy forecasts 
from every initial condition. 

The construction of a normal form transform that avoids anticipatory con- 
volutions follows the same sort of steps as described in detail in Section I2.2L 
There is no point in redoing such detail. Instead I state the coordinate 
transform as derived and checked by computer algebra that also derived and 
checked the transform leading to flHUj) and (|3"T|) . 

The stochastic coordinate transform 

x = X + ffXe- t *4)-ff 2 X(I+2e- t *)(e- t *ct)) 2 + C(e 4 ) , (41) 
y = Y + X 2 +cr[l +4Y + 2X 2 (1 - e _t ★ )] e _t ★ * 

+ o 2 [-2e~ t * +4Y(1 - 2^* )] (e^* cb) 2 + C>(e 4 ) , (42) 

transforms the system of SDEs ([5]) to the equivalent system 

X = -X 3 - XY - o-Xct) - 4aXYe _t *$ + 2o- 2 Xcbe~ t * 4> + O (e 4 ) , (43) 
Y = -(1 +2X 2 + 2Y)Y-4crY(ct> + 2Ye- t *4>) + 8a 2 Yct)e^ ^ ct> 

+ 0(e 4 ). (44) 

Note two aspects: there are no anticipatory convolutions; and X and Y vari- 
ables are different to those of Section 12.21 

The stochastic slow manifold is attractive This SDE system has Y = 
as an invariant manifold as every term in the Y SDE ( 1441) is multiplied by Y. 
Thus, in the transformed coordinates the stochastic slow manifold (ssm) is 
Y = . This SSM exponentially quickly attracts at least some finite domain 
about the origin in (X, Y, a) space as the dominant terms in the evolution 
are Y m — Y. From the stochastic coordinate transform (1421) . the SSM is 
parametrically given by (1411) and 

y =X 2 + o-[l +2X 2 (1 -e- t *)]e- t *ct>-2o-V t * [(e^**) 2 ] +£>(e 4 ) . 

On this SSM the evolution, from (1431) . is identical to the SSM model (1301) . 
This normal form coordinate transform easily displays the SSM. 
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We cannot make accurate forecasts Suppose we specify some initial 
state (Xq, Yq), either deterministic or stochastic. What forecast can we easily 
make with the SSM model ( !30l) ? In general, none. The reason is that in the 
evolution to the SSM, the SDE ( )43j) shows the slow X dynamics are coupled 
to the fast Y dynamics. But the point of deriving a slow model, for most 
purposes, is to avoid resolving the details of the fast dynamics; thus we cannot 
rationally project from (Xo, Yo) onto the SSM. In contrast, the normal form 
of a deterministic system empowers rational projecti on from nearby initia l 
conditions onto the slow model for accurate forecasts (ICox Sz Robertslll995l ). 
Abandoning Principle [2] means we cannot make accurate forecasts. 

Because it adheres to Principle [21 the stochastic normal form of Sec- 
tion \'2.'2\ similarly to the deterministic normal form, empowers rational pro- 
jection from nearby initial conditions onto the SSM. But there is a catch: in 
order to do the projection we need to anticipate the noise. Since we gener- 
ally will not know the future noise, the stochastic normal form of Section 12.21 
also cannot be used for accurate forecasting. In this sense the two stochastic 
normal forms have equivalent power. However, there is a difference. The 
anticipatory stochastic normal form of Section 12.21 has explicit convolutions 
for the projection: we may not know what they are, but we could certainly 
use the convolutions to estimate bounds and distributions for the projection 
of initial conditions. In contrast, the stochastic normal form of this section 
keeps the such information encrypted in the coupled fast and slow dynam- 
ics of 043|) and (PHI) . Consequently, maintaining Principle^ decoupling the 
slow modes from the fast, appears more powerful than avoiding anticipatory 
convolutions. 



3 Normal forms of SDEs for long term mod- 
elling 

This section uses formal arguments to establish a couple of key generic prop- 
erties of stochastic normal forms seen in the example SDE system of the 
previous section. We establish firstly that a stochastic coordinate transform 
can decouple slow modes from fast modes, to make the stochastic slow man- 
ifold (ssm) easy to see, and secondly that although anticipation of the noise 
may be necessary in the full transform no anticipation need appear on the 
SSM. 

Consider a general system of SDEs for variables x(t) £ IR m and y (t) £ M. n : 

x = Ax + f(x,y,t) , (45) 
y = By + g(x,y,t), (46) 
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where 



the spectrum of A is zero and for simplicity we assume A is upper 
triangular with all elements zero except possibly Atj for j > i (such as 
in the Jordan form appropriate for position and velocity variables of a 
mechanical system); 

for simplicity assume matrix B has been diagonalised with diagonal 
elements pi, . . ., |3 n , possibly complex, with 9^(3j < ojfl 

f and g are stochastic functions that are "nonlinear", that is, f and g 
and their gradients in x and y are all zero at the origin; 

the stochastic nature of the system of SDEs arises through the depen- 
dence upon the time t in the nonlinearity f and g — assume the time 
dependence is implicitly due to some number of independent white 
noise processes (j>k(t) (which are derivatives of independent Wiener 
processes). 



For such systems, iBoxlerl (119891 ) guarantees the existence, relevance and ap- 
proximability of a stochastic centre manifold for (HBT - H61) in some finite neigh- 
bourhood of the origin. We call this a stochastic slow manifold (ssm) because 
we assume matrix A does not have complex eigenvalues (oscillatory dynamics 
are considered briefly in Section [5]). 

For example, the toy SDE system §5§ takes the form (HoTHoT) with variables 
x = (\/o", x) and y — y , then 



A 








B 



1 





-xy 
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x 2 -2-y 2 + v^ 2 cKt). 



In principle, the matrices A and B could also depend upon the reali- 
sation of the noise. When the Lyapunov exponents of the corresponding 
linear dynamics are zero and negative respect ively, then a stochastic centre 
manifold still exists and has nice properties (iBoxlerl Il989l ). However, here 
I restrict attention to the algebraically more tractable case when the basic 
linear operators A and B are deterministic. 



Stochastic singu l ar pe rturbation systems such as those explored by 
Berglund fc Gentd (120031 ). are a subset of the systems encompassed by ( 1451 - 
^6]). For example, transform the deterministic singular perturbation system 

x = f(x,-y), y = -g(x,-y), (47) 



3 If matrix B is in Jordan form, rather than diagonalisablc, then extensions of the 
arguments lead to the same results. 
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into the form (|45ff4"S)) . First, change to the fast time T = t/e so that 



dx cry 

— = ef(x,-y), — = g[x,y). 
ax dT 

Then change to a coordinate system £, and n, where r| = is the curve 
g(x,y) =0, in which the system takes the form 

p = eF(t,Ti), ^ = P(£,)ti + G(£ ;) t 1 ). 
dT dT 

Consequently, in variables x = (v^) an d y = "H , the curve (x, y ) = (0, £,, 0) 
are a set of equilibria, at each of which the dynamics are of the form (I45T 
|4T)]) . Consequently there exists a slow manifold around each point of the 
curve, which as a whole forms a slow manifold in a neighbourhood of the 



curve (jCarrlll98ll ). Thus, the analysis presented here also applies to singular 



perturbation problems by a change in time scale and coordinate system. 



A stochastic coordinate transform We transform the SDE ( 145T - H61) in 
(x,y) to the new (X,Y) coordinate system by a stochastic, near identity, 
coordinate transform 

x = X + £,(X,Y,t) and y = Y + T|(X,Y,t) . (48) 

This stochastic coordinate transform is to be chosen such that the SDE (j4"5T 
|4T)]) transforms to a "simpler" form from which we may easily extract the 
SSM. Based upon the experience of Section l2~2l we seek to simplify the SDEs 
according to Principles p:secular-[5j and allowing anticipation. 



3.1 Transform the fast dynamics 

Suppose (H51) is some approximation to the desired coordinate transform. 
Iteratively we seek corrections and r\' to the transform, namely 

x = X+£,(X,Y,t) + £/(X,Y,t) and y = Y+rj(X, Y, t) +ti'(X, Y, t) . (49) 

Find corrections such that the corresponding updates to the evolution, say 
F' and G' in 

X = AX + F(X,Y,t) +F'(X,Y,t) , (50) 
Y = BY+G(X,Y,t) + G'(X,Y,t), (51) 

simple as possible (Principle HI). 
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For the fast dynamics, the iteration is to substitute the corrected trans- 
form (|49|) and evolution (j5Ul) - (l5*Tl) into the governing sde (|46p for the fast 
variables. Then drop products of corrections as being negligible, and ap- 
proximate coefficients of corrections by their leading order term. Then the 
equation for the jth component of the correction to the transform of the fast 
variable and the new fast dynamics is 

g ; + ^ - m; + Y_ p*y^ = x*m • ( 52 ) 

Here, Req^g| i denotes the residual of the jth component of the SDE ( 1461) . In 
constructing a coordinate transform we repeatedly solve equations of this 
form to find corrections. 

We find what sort of terms may be put into the transformation r] and what 
terms have to remain in the Y evolution by considering the possibilities for the 
right-hand side. The transform is constructed as a multivariate asymptotic 
expansion about the origin in (X, Y) space; thus all terms in the SDEs are 
correspondingly written as asymptotic expansions in (X,Y). Suppose the 
right-hand side, the residual Rea^^, has, among many others, a term of the 
multinomial form 

m n 

c(t)x^=c(t)nxrri Y r> 

t=1 j=1 

for some vectors of integer exponents p and q . Because of the special form 
of the 'homological' operator on the left-hand side of (|5*2"]) . seek contributions 
to the corrections of 

G; = a(t)X p Y q and ^ = b(t)X p Y q . 

Then this component of (l52j) becomes 

n 

a + b — \ib — c where \x — |3j — V| q £ (3 e . (53) 
Three cases arise. 

1. In the case \x = , we need to solve a + b = c where we want to put 
as much into b as possible (Principle H]). Generically, the forcing c(t) 
will have mean and stochastically fluctuating components. Neither of 
these can be integrated into b as they both give rise to secular terms 
(Principled]): the mean of c generates linear growth; the stochastically 
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fluctuating part of c almost surely generates square-root growth]^ Thus 
the generic solution is a = c and b = 0, that is, assign c(t)X p Y q to 
the Y evolution and nothing into the coordinate transform rj. 

One example when this case occurs is when all exponents in q are zero 
except for qj = 1 . Then the contribution to the Y evolution is simply 
c(t)X p Yj. Being linear in Yj, this contribution maintains Principle [3] 
that Y = is the SSM in the transformed coordinates. 

In general, since £H(3£ are all negativejf] the case \± = can only arise 
when at least one of the exponents q is positive in order for the sum 
in (1531 to be zero. Hence, generally there will be at least one Y^ factor 
in updates G' to the Y evolution, and so we maintain that Y = is 
the SSM. 

2. When < , a solution of (1531 is to place all the forcing into the SSM, 
b = e Ht * c , and do not introduce a component into the Y evolution, 
a = . As < , the convolution is over the past history of the 
noise affected forcing c(t); the convolution represents a memory of the 
forcing over a time scale of 1 /|SH(J.|. This case of *Ryi < arises when 
— 9^(3 j is large and the exponents q are relatively small, corresponding 
to low order nonlinear factors of a rapidly dissipating mode. 

However, for large enough exponents q, that is for high enough order 
nonlinear terms, the rate must eventually become positive. In the 
transition from negative to positive, the rate IHfx may become close to 
zero. Then the time scale 1 becomes large and may be as large 

as the macroscopic time scale of the slow dynamics of interest. In that 
case set the transform b = and assign this term in the forcing into 
the Y evolution with a = c . The intended use of a macroscopic model 
defines a slow time scale and consequently affects which terms appear 
in the model. 

3. When > , and accepting anticipation in the transform, we simply 
set b = e Ht * c , and do not change the Y evolution, a = . 

Consequently, we are always able to find a coordinate transform which main- 
tains that Y = is the SSM. 

4 In contrast, when the forcing c(t) is periodic, instead of stochastic, then the the 
forcing c(t) may be integrated into the coordinate transform b, instead of being assigned 
to the evolution a. 

5 Morc generally, provided the eigenvalues (3j are all non-zero whether real or complex, 
the argument still holds. Thus we can maintain Principle G2 over a very wide range of 
circumstances. 
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You may have noticed that I omitted a term in (152j) : the term g^-A^iXi 
should perhaps appear in the left-hand side. However, my omission is ac- 
ceptable when the matrix A is upper triangular as then any term introduced 
which involves X« only generates extra terms which are lower order in X^. 
Although such extra terms increase the order of X| for i > i , successive 
iterations generate new terms involving only fewer factors of X^ and so itera- 
tion steadily accounts for the introduced terms. Similarly for the Y variables 
when the linear operator B is in Jordan form due to repeated eigenvalues. 
Discussing equation (1521) for corrections is sufficient. Analogous comments 
apply to the the slow dynamics to which we now turn. 

3.2 Transform the slow dynamics 

For the slow dynamics, each iteration towards constructing a stochastic co- 
ordinate transform substitutes corrections to the transform f|4^1) and the evo- 
lution fl5*0145Tl) into the governing SDE (T4"5"]) for the slow variables. Then drop 
products of corrections as being negligible, and approximate coefficients of 
corrections by their leading order term. Then the equation for the jth com- 
ponent of the correction to the transform of the slow variable is 

3£,- n dV 

f ; + + L 3«v t ^ = Re S • ( 54 ) 

Here, Res^| i denotes the residual of the jth component of the SDE fj4"5"]) 
evaluated at the current approximation. In constructing a stochastic coordi- 
nate transform we repeatedly solve equations of this form to find corrections. 
The difference with the previous discussion of the fast variables is that the 
left-hand side of ( 1541) does not have an analogue of the — |3jT| ' term. 

We find what sort of terms may be put into the transformation correc- 
tion and what terms have to remain in the X evolution, via the correc- 
tion F , by considering the range of possibilities for the right-hand side. In 
general, the right-hand side residual Res^^ is a sum of terms of the form 

m n 

C (t)XPY"=c(t)nxrn Y f - 

for some vectors of integer exponents p and q. Because of the special form of 
the 'homological' operator on the left-hand side of f[5"4"j) . seek corresponding 
corrections 

Fj = a(t)X p Y q and fcj = b(t)X p Y q . 
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Then (1541) becomes 



a + b — \xb = c where \i = — } qifti. (55) 

e=l 

Two cases typically arisej^l 

1 . The case \± = only arises when the exponents q = as the exponents 
have to be non-negative and 91(3 f < . We need to solve a+b = c where 
we want to put as much into b as possible (Principle H]). Generically, the 
forcing c(t) will have mean and stochastically fluctuating components. 
Neither of these can be integrated into b as they both give rise to 
secular terms (Principle [T]): the mean of c generates linear growth; 
the fluctuating part of c almost surely generates square-root growth[] 
Thus at first sight the generic solution is a = c and b = , that is, 
assign c(t) YliL-\ ^-T to the X evolution and nothing into the coordinate 
transform. 

But recall Principle [5j we do not want fast time integrals in the slow 
evolution. Consider the case when the forcing has the form of a fast 
time convolution c = e wt * C . From (TT41) deduce 



■vt n 1 1 d t 1 1 dc 

c = e * C = — C H — (e * C) = — C H — . 

\y\ v dt rv| v dt 

Hence to avoid fast time memory integrals in the slow X evolution 
(Principle E]), set a = C/\v\ and b = c/v. If C(t) in turn is a fast 
time convolution, then continue the above separation. This separation 
corresponds to the integration by parts that Section 12.21 uses to avoid 
fast time, memory convolutions in the slow evolution. 

When the forcing c is a quadratic product of convolutions, then similar 
transformations and integration by parts eliminates all memory from 
the slow variables except terms of the form c(t)e vt ^c(t) where c has no 
convolutions. Algebraic transformations cannot eliminate such terms; 
for now accept the violation of Principle [5] in such quadratic forcing 
terms. 

Since the case \x — can only arises for terms in the residual with no 



6 The case £H(_l < cannot arise as all the decay rates — !ER|3j > when there are no fast 
unstable modes. 

7 But sec footnote IH pPM 
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Y dependence, we can maintain that the slow evolution of the X vari- 
ables are independent of Y, and this holds both on and off the SSmJ^I 

2. The remaining case when > occurs when at least one of the 
exponents q is positive. Accepting anticipation in the transform, we 
simply assign b = e^-kc , and do not change the X evolution, a = . 

By anticipating noise we are always able to find a coordinate transform which 
maintains a slow X evolution that is independent of whether the system is on 
or off' the SSM. Thus the projection of initial conditions and the exponential 
ap proach to a solut ion of the slow variables, called asymptotic completeness 



by iRobinsonl (119961 ). is only assured by anticipation of the noise. 

The preceding arguments are phrased in the context of an iteration scheme 
to construct the stochastic coordinate transform and the corresponding evo- 
lution. Each step in the iterative process satisfies the governing SDEs to 
higher order in the asymptotic expansions. By induction, we immediately 
deduce the following proposition. 

Proposition 1 with stochastic anticipation allowed, a near identity stochas- 
tic coordinate transformation exists to convert the stochastic system (TJJfjT^P 
into the normal form 

X ~ AX + F(X,t), (56) 
Y ~ [B + G(X,Y,t)]Y, (57) 

where ~ denotes that these are equalities to any power of (X, Y) in an asymp- 
totic expansion about the origin; there generally are exponentially small er- 
rors. 

Note: F and G may contain fast time memory integrals but these need only 
occur as products with other noise processes; for example, see (1261) and ( 12911 . 



This proposition corresponds to the general Theorem 2.1 of lArnold fc Xu Kedai 



(119931 ) except they do not identify that memory integrals may be mostly elim- 



inated. 



8 However, when the fast dynamics contain rapidly oscillating, non-decaying modes, 
the corresponding eigenvalues occur as complex conjugate pairs which typically interact 
to cause u. = ; among r apid oscillation we ca nnot completely decouple the slow modes 
from the fast oscillations ( Cox fc Roberts! fl99^ . Physically, waves do interact with mean 
flow. 
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3.3 Slow dynamics do not need to anticipate the noise 

Despite the presence of anticipatory convolutions appearing in the stochastic 
coordinate transform, we here argue that none of them appear in the slow 
dynamics because the a n ticipa tory convolutions always involve fast variables. 



Bensoussan fc Flandolil (119951 ) correspondingly show we do not need to an- 
ticipate noise on a stochastic inertial manifold. 

In the previous sections, the anticipatory convolutions only occur when 
the rate [i > . But for both the slow and the fast components, this oc- 
currence can only be generated when at least one fast Yj variable appears 
in the term under consideration. Moreover, there is no ordinary algebraic 
operation that reduces the number of Y factors in any term: potentially the 
time derivative operator might, 

d 9 v- 9 v- 9 

dt = 3t + YkB£ ' k 9Y^ + 2- XkA ^9)Q ' 

but although for non-diagonal A and B, in the algebra X« variables may be 
replaced by Xic and Y{ variables may be replaced by Y^, nonetheless the same 
number of variables are retained in each term and a Y variable is never re- 
placed by an X variable. The reason is that the x and y dynamics are linearly 
decoupled in the system (HoTfloT) . Consequently all anticipatory convolutions 
appear in terms with at least one component of the fast variables Y. 

Since the evolution (156jl of the slow modes X is free of Y variables, the 
evolution is also free of anticipatory convolutions. However, as seen in exam- 
ples, there generally are anticipatory convolutions in the evolution fl56l) of the 
fast modes Y. Further, although the stochastic coordinate transform (|48p has 
anticipatory convolutions, on the SSM Y = there are none. Consequently 
the preceding formal analysis leads to the following proposition. 

Proposition 2 although stochastic anticipation may be invoked, there need 
not be any anticipation in the dynamics ( 1561) of the slow modes in the stochas- 
tic normal form of the system ^J^4^)- Moreover, on the SSM, Y ~ ; the 
stochastic coordinate transform (EES"]) need not have anticipation. 



In contrast, lArnold fc Xu Kedail (119931 ) and lArnold fc Imkellerl (I19981 ) 



record anticipatory convolutions in the slow modes of their examples, re- 
spectively (12) and (4.6). Such anticipatory convolutions are undesirable in 
using the normal form to support macroscale models. 
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4 Implications for multiscale modelling 



This section describes some of the generic consequences of the previous sec- 
tions in modelling stochastic systems. 

Anticipation All who write down and then use coarse scale models of 
stochastic dynamics implicitly are soothsayers. In writing down a coarse 
scale model, researchers neglect the many details of any quickly decaying 
insignificant ignored modes. Proposition [1] assures us that normally this 
neglect requires us to know aspects of the near future of the ignored modes 
in order to decouple the coarse modes from the uninteresting details. In 
particular, providing initial conditions for the coarse model requires looking 
into the future. Nonetheless, Proposition [2] assures us that non-anticipative 
coarse models do exist and may be accurate for all time. 

4.1 Compare with averaging 



Papavasiliou &: Kevrekidisl (120061 ) .§5. explored the multiscale, equation free, 



modelling of the simple, two variable, one noise, stochastic system 



dx = -(y+y z )dr, (58) 

dy = — {y -x)dT+ -= dW T . (59) 
e V e 

This system has two time scales for small parameter e: for small e the fast 
variable y decays quickly to y ~ x on a t time scale 0(e); substituting 
this approximate balance into ( 158]) gives dx ~ — (x + x 2 )dT in the absence of 
noise. That is, over t times longer than 0(e) the slow variable x evolves. We 
compare the information provided by averaging to that provided by stochastic 
normal forms in multiscale modelling. 

Many apply met hods of singular perturbations to s ystems of the form 



591) . For example, iPapavasiliou fc Kevrekidisl (120061 ) use the method of av- 



eraging to deduce that 

x = x + 0( v / e) where dx = -(x + x 2 + ±)dT . (60) 

That is, solutions of (|58 l) -(l59 l) are modelled to an error 0(y/e) by the deter- 
ministic ODE (1601) which applies over T times longer than 0(e). The noise 
in the fast variable y generates the extra drift — \ dT in (IBDj) through the 
quadratic nonlinearity in the slow equation (I58p . However, averaging gives 
no basis for improving the 0(i/e) error: such errors are often large in ap- 
plications as the scale separation may only be an order of magnitude or two; 
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for example, iPapavasiliou Sz Kevrekidisl (120061 ) simulate SDEs (l58l)- f[5Tll) with 



scale separation e = 0.01 implying errors are roughly -\/e = 10% . Nor does 
averaging recognise the stochastic fluctuations, seen in the Example SDE (OQ), 
induced in the slow variable x through fluctuations in the fast variable y. 
Stochastic normal forms extract both effects, and more as well. 



As in previous sections, computer algebra modified from that by [Roberts 



(j2006d ) readily derives a stochastic normal form for the system 
But first we avoid the straight jacket of singular perturbations by simply 
rescaling time to t = r/e : that is, we adopt a time scale t where the rapid 
transients decay on a t time of C*(l), and the slow variable x evolves on 
long times At ~ 1 /e . The example system fl58|) - fl59|) is then identical to the 
Stratonovich system 

x = -e{y+y 2 ), (61) 
y = -y + x+acKt), (62) 

when the new noise magnitude c = 1 . I introduce the noise magnitude 0" 
in the SDE system ( I6ip - (I62I) in order to control truncation of noise effects. 
Minor modifications of the previous computer algebra then discovers that 
the stochastic coordinate transform 

x = X + e(Y+ ±Y 2 + 2XY) 

+ eo- [(1 + Y + 2X)e -t * cb + Ye +t * cb] 

+ ^ecr 2 (e- t *c|)) 2 + 0(e 2 ) , (63) 
y = Y + X + cre^cb-^eY 2 

+ ea[(1 - Y + e-^)e- t *(j) + Ye +t * cb] 

+ ecr 2 e- t * [cbe - ** cb + \{e~ x * cb) 2 ] +C(e 2 + X 2 ) , (64) 

maps the SDE system ( I61l) - fl6"2l into the following Stratonovich SDE system 
for the new variables X and Y: 

X = -e(X + eX + X 2 ) - ecr(1 - 2e + 2X)cb - ec^cbe^ * cb 

+ 0(e 4 + X 4 + cT 4 ) , (65) 

Y = Y[(-1 +e + e 2 + 2e 3 ) + (2e + 4e 2 )X + cr(2e + 6e 2 )cb] 

+ C(e 4 + X 4 + a 4 ) . (66) 

As before, the utility of this normal form transformation is that the SDE (166]) 
shows that the transformed fast variable Y — > exponentially quickly from 
a wide range of initial conditions for small scale separation e. Moreover, the 
methodology may refine the approximation, through further iteration, to suit 
a wide range of specified finite scale separation e. 
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This normal form transformation also shows that the new slow variable X 
evolves independently of the fast variable Y, see (KJ51) . both throughout the 
initial transient as well as thereafter: there are no initial transients in X. Fur- 
thermore, being just a transform form of the original SDE ( |62l) . the SDE (J65]) 
for the slow variable X(t) applies for all times, albeit to the truncation error; 
in contrast, the averaged system generally only applies for a finite time span. 

Although not immediately apparent, the leading approximation of the 
slow X evolution (I65jl is the averaged model (J60l) . The quadratic noise term 
in ([65]) generat es a mean drift and an effec tive new noise over long times: 



Roberts! (j2006d ) and IChao Sz Roberts! (119961 ) argued that over long times 



(j>e _t *(|)i-> 5 + 75*, (67) 

where 4>(t) is a new 'white' noise process independent of the original noise 
process cb(t). Thus the slow variable SDE ( 1651) is effectively the SDE 



X 



e{\o 2 + X + eX + X 2 ) - ea(1 - 2e + 2X)<|> - ecr 2 ^4> . (6* 



Reverting to the original (slow) time t, setting cr = 1 to match the original 
noise intensity, and re-expressing, the SDE ( 1681) becomes 



dX = -(I + X + X 2 + eX)dT- v/e(1 -2e + 2X)dW T - J\ dW T . (69) 



The deterministic ode found at leading order, dX = — [X + X + X 2 )dT, is 
the averaged model ( l60l) . However, the SDE ( l69|) also makes explicit some 
of the errors in averaging. Firstly, the y/e error of the averaged model ( 1601) 
comes from its neglect of the stochastic fluctuations: to leading order we can 
combine noise processes YV and W to determine that the slow variables are 
better modelled by the SDE dX = -{\ +X + X 2 )dT+ y/3e/2dW T ; although 
the two stochastic terms in the SDE (IBT)]) are even better. Secondly, the 
SDE (1691) also discerns 0(e) contributions to the deterministic terms which 
may well have significant effects at finite scale separation e. Simple averaging 
misses all of these effects. 



4.2 Equation free simulation 



Kevrekidis et al.l (120031 ) promote a framework for computer aided, equation 
free, multiscale analysis, which empowers systems specified at a microscopic 
level of description to perform modeling tasks at a macroscopic, systems 
level. When the microscopic simulator is stochastic, that is Monte-Carlo, 
or effectively stochastic, such as molecular and discrete element simulators, 
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then the issues addressed in this article of the nature and extraction of long 
term dynamics from a stochastic system are crucial to the equation free 
methodology 

Equation free modelling is designed to solve specific multiscale systems 
with specific finite scale separations. Thus a challenge for future research is 
to maintain reasonable accuracy in estimating the long term dynamics by 
extracting from numerical realisations the sort of information extracted by 
these normal form coordinate transforms and without knowing any algebraic 
representations of the systems of interest. The stochastic normal form trans- 
formation shows what might be achieved in principle. The challenge is to find 
out how to achieve it from just a finite number of short bursts of realisations. 

On the macroscale the stochastic effects may be relatively small. How- 
ever, a deterministic macroscale model is often structurally unstable: one 
example is the structural instability of the averaged model dx = dt for the 
SDE (JTJ); instead we prefer the stochastic model dX = e dW of the SDE ([31). 
Moreover, even on the macroscale a deterministic model for some averaged 
slow variable is almost inevitably different from the average of the system 
with noise included. This difference follows from the same line of argu- 
ment that establishes that the expectation of realisations is generally dif- 
ferent from the expected position of the stocha stic slow manifold, see fl34l ). 
Noise induced mean drift must be r ecognised (ISri Namachchivaya fc Leng 
199(1 ISri Namachchivava fc Linlll99lh . 



To do coarse step integration we need to estimate the macroscale drift 
from the microscale simulations. Noise hides this drift making accurate esti- 
mation difficult. We need either long bursts of microsimulations, or many re - 
alisations, or variance reduction techniques (IPapavasiliou fc Kevrekidisll2006l ). 
or a combination of all three. 

As well as the drift, the fluctuations in the macroscopic quantities should 
be modelled. Thus the macroscale integration should be that of a system of 
SDEs. Because the macroscale SDEs model microscale processes, I conjecture 
that the macroscale SDEs must be interpretted as Stratonovich SDEs. The 
challenge is to develop Stratonovich integrators that only use short bursts of 
realisations. 

In equation free simulation one projects a macroscopic time step into the 
future, then executes a burst of microscale simulatio n in order to estimate the 
macroscopic rate of change (IKevrekidis et al.l 120031 ). Initial rapid transients 
must be ignored in each burst as the microscopic system attains a quasi- 
equilibrium. In a stochastic system, the true SSM can only be identified via 
integrals over fast time scales, see Section 12.41 However, these are generally 
integrals of both the past and the future. Thus, to estimate macroscopic rates 
of change in a stochastic system, we must not only neglect initial transients, 
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but also data from the end of a burst of microscopic simulation in order to 
be able to account for the integrals which anticipate the noise processes. 

Lastly, the gap-tooth scheme empowers equat ion free modelling across 
space scales as well as time scales ( iGear et al.ll2003L e.g.). For spatiotemporal 
stochastic systems we need theoretical support for the notion that spdes can 
be modelled by the gap-tooth sch eme in the same way as deterministic pdes 
([Roberts fc Kevrekididlto appear! ) . Only then will we be assured that we can 
cross space scales as well as time scales. 



5 Long time modelling of stochastic oscilla- 
tions 

The previous sections focus on the separation of slow modes from fast modes 
in stochastic systems. Persistent oscillations are another vitally important 
class of dynamics. Hopf bifurcation is the example considered in this section, 
but many other cases occur including wave propagation. The challenge ad- 
dressed here is how to consistently model the evolution of oscillations over 
long time scales when the oscillations are fast and in the presence of stochas- 
tic noise fluctuations over all time scales. To model over long time scales we 
eliminate from the model all fast time dynamics. 

As an exampl e let us explore the stochast ic P uffing-van der Pol dynami cs 
also analysed by lArnold fc Xu Kedail (119931 ) and lArnold fc Imkellerl (119981 ): 



X! = [oc + crcb(t))x-| + (3x) — x^ 



X^X] 



(70) 



where , as before, (J) is some white noise process. lArnold. Sri Namachchivaya fc Schenk-Hoppe 
(119961 ) describe the importance of the stochastic system (1701) in applications. 
In the absence of noise, cr = , this system exhibits 

1. a deterministic pitchfork bifurcation as the parameter <x crosses zero 
with [3 fixed, say |3 = — 1 for defmiteness; and 

2. a deterministic Hopf bifurcation as the parameter |3 crosses zero with 
<x fixed, say <x = — 1 for defmiteness. 



In the presence of noise, cr > 0, computer algebra ([Roberts! 12006 a ) readily 
derives the stochastic normal form for the Duffing-van der Pol equation ( 1701 
near the stochastic pitchfork bifurcation when parameter <x crosses zero with 
fixed (3. This section explores the issues arising when constructing a normal 
form for the Duffing-van der Pol equation (I70p near the stochastic Hopf 
bifurcation as the parameter (3 crosses zero with cx = — 1 . 
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Figure 4: stochastic Hopf bifurcation in the Duffing-van der Pol oscilla- 
tor (ITUl) as parameter (3 crosses zero with a = — 1 and noise amplitude 
a = 0.5 . Two realisations are plotted in each case. 



5.1 Approaches to stochastic Hopf bifurcation 



Figure H] shows the noisy Duffing-van der Pol oscillator (1701) of lArnold fc Imkeller 
(119981 ) for parameter (3 = ±0.1 . Deterministically (cr = 0), as parameter (3 
crosses zero, a Hopf bifurcation occurs from the stable equilibrium at the 
origin to a stable limit cycle with frequency 1 . What happens in the pres- 
ence of parametric stochastic forcing when <y ^ 0? Figu re H] reaffirms that 



a nois y version of the stochastic bifurcation takes place. ICoullet fc Spiegel 
(119831 ) first explored a normal form of Hopf bifurcations with noise. Fur- 
ther research on such sto chastic bifurca t ions e lucidated some fascinating fine 
structure. For example, Keller fc Ochsl (119991 ) explored the structure of the 
random 'lim it cycle' attractor using a s tochastic yersion of the subd i vision 
algorithm of iDellnitz fc Hohmannl (119971 ); whereas lArnold fc Imkellerl ( 119981 ) 
explored the structures using a normal form approach very close to that used 
here. However, the emphasis here is not on the stochastic Hopf bifurcation 
as such, but instead using it as the simplest prototype system with stochastic 
oscillatory dynamics. We look at the issues afresh to explore the characteris- 
tics of a long term stochastic model of such stochastic oscillatory dynamics. 
In the future, these considerations will underpin the multiscale modelling of 
stochastic oscillations and waves. 

Section 15.21 constructs a stochastic coordinate transform from which we 
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may extract significant properties of a stochastic Hopf bifurcation. Solutions 
of the Duffing-van der Pol oscillator ([TO]) are most conveniently represented 
in complex exponentials as 

xt « a(t)e it + b(t)e- it , (71) 

where for real solutions Xi, the amplitudes a and b are complex conjugates. 
Then Section 15.21 finds a stochastically forced Landau model governing the 
evolution of the complex amplitudes a and b: 

d « |(3a-(^-|i)a 2 b + o-v^72(a4)o-bct) +2 ), (72) 
b « ipb-(| + |t)ab 2 +aV672(b(l)o-ac|)_2), (73) 

to errors (9((3 2 + a 2 + e 4 ) where e = |(a, b)| measures the size of the oscilla- 
tions, and where cb m (t) are independent 'white' noises arising from the forc- 
ing components near frequencies and ±2 in the applied noise process 4>(t); 
'near' means within ±5 of the specified frequency. This model resolves the 
slow evolution of the complex amplitudes near the Hopf bifurcation, small (3, 
under the influence of the nonlinearity and a weak stochastic forcing, small cr. 
This model empowers long term simulations with efficient large time steps 
as the complex amplitudes are slowly-varying. 

Note: the analysis also applies in the case of the deterministic forcing 4> = 
cos2t, for which <t)o = and \/28<$)±2 = \. Then the above model, d ~ 
^(3a+^crb and its conjugate, successfully predicts the Mathieu-like instability 
with eigenvalues A = ^(3 ± \a . 

View ( ITT]) as a time dependent coordinate transform of the (xi , X] ) phase 
plane. In principle, any dynamics in the phase plane may be described by 
the evolution of the complex amplitudes a and b. The utility of the coordi- 
nate transform ( ITTI) is that it empowers a simple description of oscillations 
with frequency near 1 : namely (l72l)-(173l) for the Duffing-van der Pol oscilla- 
tor ( ITOl ). However, to simply describe such nonlinear stochastic oscillations 
Section 15.21 modifies the coordinate transform (ITT]) through nonlinear and 
stochastic terms. That is, there is a time dependent, coordinate transform 
of the phase plane that leads to the normal form (|T2]) -(IT3 l . 

I emphasise this different view of ( ITT]) . Many would view (ITT]) as an 
approximation to x(t) that can only resolve slowly varying oscillations. In 
contrast, I present (ITT]) as the leading term in a coordinate transform, a 
reparametrisation, of the entire phase (xi,xi) plane that in principle en- 
compasses all dynamics in the phase plane. The approximate model then 
arises by finding parameter regimes, in this new coordinate system, where 
the evolution of 'coordinates' a and b is usefully slow. 
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Amplitude/phase models do not decouple lArnold fc Imkellerl (119981 ) 
analysed a Hopf bifurcation by transforming to real amplitude r and phase 
angle (p coordinates and deducing a model f = • • • and (p = 1 + • • • . This 



appro ach is certainly effective for unforced deterministic problems ([Roberts 
2006a . e.g.). However, the presence of time dependent forcing, whether 
stochastic or deterministi c, breaks time translation symmetry. Consequently, 
Arnold &: Imkellerl (119981 ) must couple the phase c p back into th e ampl itude r 
evolution, as also seen in the normal form (39) of I Arnold et al.l (119961 ). Such 
coupling of the fast phase into the notionally slow amplitude confounds our 
aim to use the normal form for lo ng time mode l ling. 

Because of their different aim, lArnold et al.l (119961 ) convert back to a pair 
of fast Cartesian variables to obtain a canonical system that is generic for the 
class of stochastic Hopf bifurcations; thus they establish that the pattern of 
behaviour they explore is generic for Hopf bifurcations. But our aim here is 
different: we aim to construct models suitable for exploring long time evolu- 
tion; our normal form is consequently different. We use complex amplitude 
coordinates, the a and b seen in (1721) and (173"]) . as originally proposed by 
Coullet k Spiegell (Il983l ). 

Stochastic averaging seems to suffer the s ame defect of not recognising the 



broken time symmetry (lArnold et al.lll996l . equations (16-20)). Stochastic 



averaging also does not appear to detect the split in Lyapunov exponents 
present in stochastic Hopf bifurcations. 



Prefer a strong model lOlarrea Sz de la Rubial (119961 ) comment that "When 



the reduction to the normal form is done . . . only the deterministic part of 
the equations retain the characteristic radial symmetry." and then assert 
"This makes it necessary to work with the two-dimensional probability dis- 
tribution." Thus they introduce early in their analysis some probability dis- 
tributions governed by Fokker-Planck equations and hence derive only weak 
models. In contrast, here we maintain strong modelling of each realisation 
of the noise. We avoid weak models. 



5.2 Construct a stochastic normal form 

To construct the stochastic normal form for the Duffing-van der Pol oscil- 
lator (I70|) . with parameter a = — 1 , I use an iterative scheme to construct 
a useful nonlinear coordinate transform. The coordinate transform must be 
time dependent to adapt to both the oscillations and to the stochastic ef- 
fects. The starting approximation to the linear time dependent coordinate 
transform is (1711) . Iterative modifications to (ITTj) result in a description of 
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the Duffing-van der Pol oscillator (170|) which only has slow processes suitable 
for long time simulation. 



The homological equation Each step in the iteration improves the nor- 
mal form description of the dynamics. Suppose that at some step in the 
iteration, the coordinate transform and consequent evolution is 

xi=£,(a, b,t) where d = g(a, b,t) and b = h(a,b,t) 

for some known functions £,, g and h. Seek small corrections, denoted by 
dashes, to £,, g and h so that 

xt = £, + £/( a, b, t) where d = g + g'(a, b, t) and b = h + H'(a, b, t) 

(74) 

better satisfies the Duffing-van der Pol oscillator (170 p . We measure how 
well the Duffing-van der Pol oscillator (ITU]) is satisfied by its residual, Rea^jj. 
Substitute ( 1741 into the Duffing-van der Pol oscillator ( 1701) . omit products 
of small corrections, approximate £, ~ ae xt + be~ xt and g ~ h. « (3 ~ 
a m whenever multiplied by a correction, and deduce that in the complex 
amplitude coordinates, the homological equation is 

£4 + £,' + (i2g' + g ;)e lt + (-i2h/ + h^e'* + Re^ = . 

But there is one further refinement: we aim for d = g and b = b to only 
possess slow dynamics; thus, presuming this aim is possible, also omit the 
time derivatives g' t and h\( to give the homological equation 

£4 + £/ + i2g'e tt - HWe^ + Resjn] = • (75) 

This approach avoids difficulties that appear in the homological equation for 
amplitude-phase models. The homological equation ( 1751) governs corrections 
to the complex coordinate transform. 



5.3 Linear noise effects 

An iterative scheme to find a stochastic coordinate transform and c orrespond- 
ing evolution was coded into computer algebra ( jRobertsl 12006 d ). Iterative 
improvements to the coordinate transform and the model continue until the 
residual of the Duffing-van der Pol oscillator (1701) reaches a specified order 
of error. To effects linear in the noise magnitude a the iteration finds the 
stochastic model ( 1721) and (1731) to the specified errors. In terms of the Fourier 
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transform (b(Q) of the noise, cb(t) = J" e l0t cb(n) dQ, the corresponding 
stochastic complex coordinate transform is 

x^ = ae it + be- it + +i)a 3 e i3t +(1 -i)b 3 e - i3t ] 

-e i(a+1)t *(0)dO 



o~ia 



+ ffib 
+ y/25a 



D a(a + 2) 
i 



Ma-m 



d O(Q-2) 



ct>(n) da 



i(ach) - b^je 1 - ^(bcb - acj)_ 2 )e 



—it 



- |acb 2 e i3t + ibc[)_ 2 e- i3t ] + G(|3 2 + a 2 + e 4 , 5 3/2 ) , (76) 

where the integration domain D avoids singularities in the integrand as ex- 
plained in Section 15.3.21 



5.3.1 Deterministic effects 

The first line of (1761) describes the well established deterministic shape of the 
limit cycle in the deterministic Hopf bifurcation. When the residual ReqyQ] 
has terms with factors e imt for some integer m, |m| ^ 1 , and no other explicit 
time dependence, then as usual we update the complex coordinate transform 
by a correction £,' proportional to e 1TrLt /(m 2 — 1), and do not change the 
evolution, g' = h' = . 

Deterministic terms in the residual with factors e ±vt , and no other explicit 
time dependence, such as the term (i(3a + 3a 2 b)e u , are resonant and as usual 
must be assigned to correct the evolution, upon dividing by the ±2i factor of 
the homological equation (1751) ; see the deterministic nonlinear and (3 terms 
in the model fl72l)-(lT3l. 



5.3.2 Non-resonant fluctuations 

The second two lines of the transform (1761) describe how stochastic fluctu- 
ations non-resonantly perturb the oscillating dynamics. These arise from 
terms in the residual Res^Q] of the form 



cKt)e 



e i(a±1)t cMn)dQ 



Away from resonance, namely in the domain D = IAU me {_2,o,2}[ Tn— ^> m +S] , 
these terms in the residual generate the desingularised integrals in ( 1761) rl 

9 In analyses to higher order in the oscillation amplitude more resonant frequencies 
occur; for example, integrals arise with singularities at frequencies Q. = ±4 in some terms 
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Figure 5: schematic plot of three realisations of the amplitude of resonant 
noise cbo(t) for some mesoscale cutoff 5. 



Rewriting these integrals as a convolution f(t) *(j)(t) recognise that formally 
f = e ±vt J* D o(£5±2) eiat dXi . This integral for the convolution kernel f may 



be written in terms of the Sine integral (lAbramowitz fc Stegunlll965l . §5.2) 
from which we deduce that the convolution kernel f(t) decays like 1/(6|t|) 
for large |t|. Assuming that convolutions of f(t) with stochastic white noise 
do converge in some sense, the complex transform appears to necessarily 
involves the entire past and future of the noise. In contrast to the pitchfork 
bifurcation, which only needs to look a little way into the future and the past, 
in the Hopf bifurcation we look far into the future and the past in order to 
construct the stochastic coordinate transform. 



In contrast, ICoullet fc Spiegel! (119831 ). in their equations (18) and (19) 



assign the entire integral to the evolution ( 1721) - (1731) . rather than to the trans- 
formation, just because one frequency is resonant. This approach seems in- 
consistent in the neglect of the time derivatives g' t and h' t in the homological 
equation as such derivatives are large for 'white' noise. Their assignment to 
the evolution is consistent when the noise cb(t) has a narrow band spectrum 
around the resonant frequencies. 



of 0(e 2 o~). In such higher order analyses the domain of integration D will have further 
intervals excised to avoid resonances. 
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5.3.3 Resonant fluctuations 



The excised parts of the integrals in the transform ( 1761) correspond to res- 
onances. These resonances generate terms in the model (l72 i) -(173 j) involving 
components of the (complex) noise process 



4>m(t) 



i rm+5 

e i(n - m)t $(a)da, (77) 

m— 6 



26 J 



normalised so that E [|(j) m | 2 ] = 1 under the original white noise assumption 

that E [4>(Q] 4>(-Q)] = 6(H — CI) (here 5() denotes the Dirac delta function 
and * the complex conjugate); Figure [5] plots three realisations. Being a 
narrow band integral (with the dominant frequency accounted for by the 
g-imt f ac t or ) the 4>m.(t) are slowly varying noise processes: Figure shows 
cbo(t) has slow variations on the fast times scale At = 2n of the oscillations. 
They are independent of each other as the domains of integration do not 
overlap (for small cutoff 5). Each <t> m (t) has autocorrelations which decay 
on a time scale of order 1/6, roughly the width of the window in Figure 0, 
but for time scales 3> 1 /6 the autocorrelation is zero and the cb m look like 
white noise processes. Thus choose the 'cutoff' 1/6 to be a mesoscopic time 
scale: one longer than the period of the limit cycle; but much shorter than 
the long macroscopic time scale on which the model (JT2l) -( JT3l) is to be used. 
Then cb m (t) are effectively independent white noise processes in the long 
term model. 

Encouragingly, although the Fourier transform cb(O) requires the entire 
history of the noise, the parts of it that appear in the model (}T2i) -( jT3"j) are 
essentially local in time. That is, as for non oscillatory dynamics, the long 
term model itself does not require anticipation of the noise. 

The fourth and fifth lines in the transform (1761) arise through the ex- 
cision of the resonant parts of the frequency domain from the integrals in 
the coordinate transform ( 1761) . These resonant frequencies not only affect 
the evolution but also the coordinate transform as seen in these two lines 
of (J76J). 

These resonant fluctuations also force the complex amplitudes a and b 
to change their meaning in the presence of noise. I do not precisely and 
explicitly define the complex amplitudes a and b; implicitly they are the 
component in e ±xt in the oscillations. However, whatever definition one may 
try to adopt, implicitly or explicitly, the noise changes the definition through 
the terms appearing on the fourth and fifth lines in the transform (1761) . Recall 
that in non-oscillatory systems noise also changes the presumed definition of 
slow variables: for two examples, the SSMs (j2J) and ( 1321) show that we can- 
not parametrise a SSM in terms of the original slow variable x, but a new 
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variable X which is necessarily different in the presence of noise. Similarly 
here: in the presence of noise, the coefficient of e lt in the stochastic coordi- 
nate transform is not just the complex amplitude a but instead is approxi- 
mately a + \<j\flb | ( acj)o — bcp2) , and analogously for the coefficient of e~ xt . 
Noise affects the meaning of the complex amplitudes. 

These terms of the fourth and fifth lines in the transform (!7H|) . and the 
corresponding terms in the model (I72I) - (I73I) . are proportional to \fb where 5 is 
the small width of the domain excised from frequency space about the reso- 
nant terms. Can these terms be ignored as small? I contend it depends upon 
the use of the slow model (1721) -(1731 ) . In a long term simulation we may use 
macroscopic time steps of size At, say, in numerically integrating (17121 - (1731 . 
In this numerical integration we would treat the cb m (t) noises as white; thus 
their decorrelation time 1 /5 must be less than the numerical time step At. 
That is, a lower bound for the excised mesoscale cutoff is 6 > 1 /At. Thus, 
a stochastic time integrator could treat these terms as of OiA /-s/At) but no 
smaller. 



5.4 Quadratic noise effects 



In many applic ations, quadratic noise effects gene r ate important mean deter- 
minis tic drifts (ISri Namachchivaya fc Lend Il990l . ISri Namachchivaya fc Lin 
199ll . e.g.). This is easily seen in some exam ples, even using the method 
of averaging (jPapavasiliou fc Kevrekidisl 12006k §5, e.g.). Such mean drifts 



are often important. Thus, generically we must also explore how to analyse 
quadratic noise effects. 



5.4.1 Double integrals of noise complicate 

For oscillatory dynamics, as in the Hopf bifurcation of the Duffing-van der 
Pol oscillator (1701 . the outstanding complication is the appearance of dou- 
ble integrals across all frequencies in the stochastic fluctuations. Quadratic 
noise effects not involving such double integrals are straightforwardly han- 
dled as before. Terms of 0(o 2 ) will contain double integrals of the form 
Jd Jd where both O and £1 represent noise frequencies. The (black) 

hatched region in Figure E] shows this domain of integration. However, in 
the Hopf bifurcation of the Duffing-van der Pol oscillator (170 p . the kernel of 
such double integrals also has a singularity along the line O + D. = . Thus 
excise the (blue) diagonal strip shown in Figure [6] to remove the singularity 
to leave an integral over non-resonant effects in the domain D. Then addi- 
tionally analyse the excised strip as a resonant effect that directly influence 
the evolution of complex amplitudes a and b. 
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a 




Figure 6: the integration domain D x D, hatched, also has a further resonant 
region, the diagonal blue strip, excised to give the integration domain D for 
double integrals over the noise frequency. 
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Recall we use the residual of an SDE system to drive corrections to the 
normal form stochastic coordinate transform. In the residual of the Duffing- 
van der Pol oscillator ( ffOl) quadratic noise terms arise of the form 



D 



e i(n+A±1)t K±(a,a)$(a)$(a) da da 



D 



where the integrand kernels are 



{Q. + Q.±Q.Q.)[C1 + a ±2) 

2fa±2ifa±2iaa 



(78) 



Before excising the blue strip in Figure [6] to avoid the division by zero near 
a + a = , change the parametrisation of the integration domain to cu = 
— a) and cu = a + a so that O. = cu + j<v , D. = — cu + \<x> , and the 
Jacobian of the transform is one: parameter cu measures the distance from 
resonance. In this new parametrisation, the integration kernels 



K 



2(4cu 2 ^4cu - cu 2 )(2± cu] 



± 



(2cu ±4 + cu)(2cu ^4 
1 



(cu + 2)(cu-2) 



as cu 



cu)(2cu + cu)(2cu 
-> 0. 



cu 



(79) 



Then the double integrals in the residual are split into non-resonant and 
resonant parts: 



I- 



e i(n+ ^ 1)t K ± (a,a)$(a)^(a) da da 



D 



r-6 



+ e 



e lu, \Mcu) dcu, 



(80) 



-6 



where 



\i>±(cu) 



D 



and where domain D = D x D without the resonant strip as excised in 
Figure El The non-resonant double integral in flHUl) contributes components 
to the stochastic coordinate transform. The resonant integral on the second 
line of ( IHUl) contributes a component to the evolution in the new coordinates. 
Although the details will differ, the above integrals will appear in the analysis 
of general stochastic Hopf bifurcations. 

The stochastic dynamics in the normal form coordinates will involve the 
integral (l8Tl) . The integral ( IHTl) specifies the Fourier transforms of two com- 
plex conjugate components ipidO that express a nonlinear combination of 
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Figure 7: one realisation of the complex quadratically generated 'noise' 
"4>±(t) « 0.87ib r (t) ± i0.20\pi(t) where the real part is the larger blue curve 
and the imaginary part is the smaller green curve. The resonant window size 
6 = 0.2. 

the original noise process 4>(t). Here write these in terms of the real and 
imaginary parts 

*Mt) = c T ^ T (t)±ic*lH(t), (82) 

where the constants c T and Ci are chosen so the variances E [ib^:] = E Pip?] = 
1 ; these constants do not seem to vary significantly with mesoscale cutoff 6. 
Figure [7] shows one realisation of i]j±(t) illustrating that ip±(t) vary slowly 
over one period of the microscale limit cycle, and that they look like white 
noise processes over the long time scales resolved by the complex amplitudes 
a and b. In the Hopf bifurcation of the Duffing- van der Pol oscillator (170]) 
the processes ib±(t) appear to have zero mean; this may not hold for other 
stochastic Hopf bifurcations. 



5.4.2 Refine the normal form transformation 



Separa ting the double integrals as described, computer algebra ((Roberts 



2006d ) iteratively refines the stochastic coordinate transform ( 176]) and si- 
multaneously derives the following SDEs for the evolution of the complex 
amplitudes of the Duffing-van der Pol oscillator (1701) : 

d « ±|3a-(±-fl)a 2 b + av/^(acbo-b4) +2 ) 
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+ i±ff 2 (c r \b T + ic i \h)a-iSff 2 (±cb 2 1 + J(M>-2)a, (83) 
b « ^(3b-(^ + |i)ab 2 +av^72(bct)o-act)_ 2 ) 

- i±o- 2 (c r ib T - ic^Ob + iSa^cb 2 , + |cb 2 cb 2 )b . (84) 

The order of error in these SDEs is 0(e 4 + cr 3 + (3 2 , 5 3 ^ 2 ) . These SDEs account 
for more noise interactions than the lower order model (]72p - (|73p and thus 
are more accurate. 

For very small mesoscale cutoff 5, that is for simulations on very long 
time scales, the quadratic noise effects involving ib r and i[h are the dom- 
inant influences on the complex amplitudes a and b of the oscillations of 
the Duffmg-van der Pol oscillator ([TO]) . These two noise processes, see the 
integral ( |8~lT) . arise as integrals of quadratic terms in the original noise pro- 
cess cf>. Analogo usly to the qua dratic noise processes analysed on stochastic 
slow manifolds ( Roberts 2006(1 §5), as used in (1671) . I conjecture that ij> 



r- 



il>t and (t) are effectively independent when sampled over long time scales. 
Consequently, over very long time scales, one would model real dynamics of 
the stochastic Duffmg-van der Pol oscillator ( fTDl by the Stratonovich sde 



da « [1 3 a - {{ - |i)|a| 2 a] dt + i^c r a 2 a dW T - \ci0 2 a dW t , (85) 

where complex a measures the amplitude and phase of the oscillations, 
W r and Wt denote independent Wiener processes, and c r ~ .87 and Ct ~ .20 
(Figure [TJ). 

For medium mesoscale cutoff 6 use the more complete SDE model ( |83l) . 
This model, with its effects in \fb and 6, will be needed when the desired 
time resolution of a numerical simulator, essentially the integrator's time 
step At, is within a few orders of magnitude of the natural period of oscil- 
lations, here the period is about In. A challenge for future research is to 
construct special SDE numerical iteration schemes when, as here, the sde 
itself depends upon the chosen time step At; I am only a ware of SDE schemes 



which assume the SDE is indep endent of the time step ( iHigham fc Kloeden 



20051 . iKloeden fc Platen! 119921 . e.g.). Physically, the dependence upon the 
macroscopic time step is due to the difficulty in discerning what is and is not 
a resonant forcing of the oscillations, see Sections I5.3.2f[5~3~3l In multiscale 
modelling, as shown here, the macroscopic system, whethe r expressed as al 



gebra ic equations or solved using equation free methods (IKevrekidis et al. 



20031 1. may depend upon the the length or time scale chosen for simulation. 

The specific equations and formulae in the section are specific to the 
Duffmg-van der Pol equation (ITU]) . Nonetheless, I contend that the nonlinear 
and stochastic nature of these Dufnng-van der Pol oscillations are generic 
for most of the interesting issues discussed in this section. Consequently, I 



45 



conjecture that almost all long time scale modelling of stochastic oscillations 
has to address and resolve the issues discussed in this section. 



6 Conclusion 

Stochastic coordinate transforms illuminate the modelling of multiscale stochas- 
tic systems. Being a coordinate transform, a resultant 'stochastic normal 
form' describes the complete dynamics of the original system, Proposition [U 
From the normal form we easily extract the stochastic slow dynamics that 
are of interest over macroscopic times, from the uninteresting fast dynamics 



(jArnoldll2003l . §8.4, e.g.). This approach is more powerful than averaging as 



the coordinate transform may be system atically refined, especially with the 



aid of computer algebra (jRobertsI l2006d ) , and so errors are more controlled 



In contrast to earlier work, this article argues that two modelling simpli- 
fications may always be achieved without sacrificing fidelity with the original 
stochastic system. Firstly, the stochastic slow manifold and the evolution 
thereon need not have any terms anticipating the original noise processes, 
Proposition El Secondly, effects linear in the noise processes in the evolu- 
tion on the stochastic slow manifold need not involve any memory integrals 
either, Proposition [TJ Section [2] illustrates these principles for the example 
SDE system flU). 

A challenge for future research is to let the algebraic techniques used 
herein inspire development of numerical techniques useful for multiscale com- 
putations. From a finite number of bursts of stochastic realisations we need 
to determine information to empower making macroscale time steps while 
remaining faithful to the underlying stochastic dynamics. 

Section [3] explored oscillatory dynamics in the stochastic Duffing-van der 
Pol equation (1701) . It demonstrates that transforming the SDE to a slow 
model for the complex amplitude is a delicate process that requires careful 
treatment of noise integrals in order to form a consistent model of the long 
term evolution. The specific and formal analysis herein needs to be extended 
to generic oscillatory systems to discover general modelling principles. 
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